1 Setup

1.1 Import necessary R packages

library(ComplexHeatmap)
library(SingleCellExperiment)
library(circlize)
library(tidyverse)
library(DT)

1.2 Load scRNA-seq datasets:

pbmc_lupus <- readRDS("~/popuDE/reproducibility/batch_correction/pbmc_lupus_bc_dr.rds")
pbmc_covid <- readRDS("~/popuDE/reproducibility/batch_correction/pbmc_covid_bc_dr.rds")

2 PBMC-SLE

unique.celltypes <- unique(colData(pbmc_lupus)[['celltype']])
unique.subjects <- unique(colData(pbmc_lupus)[['subject']])
tab <- matrix(NA, nrow = length(unique.subjects), ncol = length(unique.celltypes),
              dimnames = list(unique.subjects, unique.celltypes))

for (celltype in unique.celltypes){
  for (sub in unique.subjects){
  tab[sub, celltype] <- sum(colData(pbmc_lupus)[['subject']] == sub & colData(pbmc_lupus)[['celltype']] == celltype)
}
}



top_barplot <- columnAnnotation(celltype = anno_barplot(x=colSums(tab),
                                                        add_numbers=T,
                                                        numbers_rot=0,),
                                height = unit(3, "cm"))
right_barplot <- rowAnnotation(subject = anno_barplot(x=rowSums(tab),
                                                      add_numbers=T,),
                               width = unit(3, "cm"))


cell_func <- function(j, i, x, y, width, height, fill){
  #if(tab[i, j] > 100)
    grid.text(sprintf("%.0f", tab[i, j]), x, y, gp = gpar(fontsize = 10))
}


Heatmap(tab, name = "Number of cells", 
        col = colorRamp2(c(0, 200, 1500),c("green","white","red")),#viridis(100), 
        rect_gp = gpar(col = "white", lwd = 2),
        column_names_rot = 30, cell_fun = cell_func,
        show_column_dend = FALSE, top_annotation = top_barplot,
        show_row_dend = FALSE, right_annotation = right_barplot,
        width = unit(9, "cm"), height = unit(20, "cm"))

composition.tab <- tab / rowSums(tab)
composition.tab %>% 
  round(3) %>% 
  datatable(extensions = 'Buttons',
            caption = "PBMC-SLE",
            options = list(dom = 'Blfrtip',
                           buttons = list('copy', 'pdf', 'print',
                                          list(extend = 'csv', filename = "PBMC-SLE"), 
                                          list(extend = 'excel', filename = "PBMC-SLE")),
                           lengthMenu = list(c(10,25,50,-1), c(10,25,50,"All"))))

3 PBMC-COVID19

unique.celltypes <- unique(colData(pbmc_covid)[['majorType']])
unique.subjects <- sort(unique(colData(pbmc_covid)[['sampleID']]))
tab <- matrix(NA, nrow = length(unique.subjects), ncol = length(unique.celltypes),
              dimnames = list(unique.subjects, unique.celltypes))

for (celltype in unique.celltypes){
  for (sub in unique.subjects){
  tab[sub, celltype] <- sum(colData(pbmc_covid)[['sampleID']] == sub & colData(pbmc_covid)[['majorType']] == celltype)
}
}

top_barplot <- columnAnnotation(celltype = anno_barplot(x=colSums(tab),
                                                        add_numbers=T,
                                                        numbers_rot=0,
                                                        height = unit(3, "cm"))
                                )
right_barplot <- rowAnnotation(severity = anno_block(gp = gpar(fill = 5:6), 
                                                     labels = c('control', 'severe/critical'),
                                                     width = unit(1, "cm")),
                               empty = anno_empty(border = FALSE,
                                                  width = unit(0.1, "cm")),
                               sample = anno_barplot(x=rowSums(tab),
                                                     add_numbers=T,
                                                     width = unit(3, "cm"))
                               
                               )


cell_func <- function(j, i, x, y, width, height, fill){
  #if(tab[i, j] > 100)
    grid.text(sprintf("%.0f", tab[i, j]), x, y, gp = gpar(fontsize = 10))
}

row_split <- colData(pbmc_covid)[c('sampleID', 'characteristics: CoVID-19 severity')] %>%
  as_tibble() %>% distinct() %>% arrange(sampleID) %>% 
  column_to_rownames('sampleID')


Heatmap(tab, name = "Number of cells", 
        col = colorRamp2(c(0, 1000, 6000),c("green","white","red")),#viridis(100), 
        rect_gp = gpar(col = "white", lwd = 2),
        column_names_rot = 0, column_names_centered =T, 
        cell_fun = cell_func,
        show_column_dend = FALSE, top_annotation = top_barplot,
        show_row_dend = FALSE, right_annotation = right_barplot, 
        row_split = row_split,
        width = unit(12, "cm"), height = unit(20, "cm"))

composition.tab <- tab / rowSums(tab)
composition.tab %>% 
  round(3) %>% 
  datatable(extensions = 'Buttons',
            caption = "PBMC-COVID19",
            options = list(dom = 'Blfrtip',
                           buttons = list('copy', 'pdf', 'print',
                                          list(extend = 'csv', filename = "PBMC-COVID19"), 
                                          list(extend = 'excel', filename = "PBMC-COVID19")),
                           lengthMenu = list(c(10,25,50,-1), c(10,25,50,"All"))))